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ABSTRACT 

The Wilkinson Microwave Anisotropy Probe (WMAP) team subtracted 
template-based foreground models to produce foreground-reduced maps, and 
masked point sources and uncertain sky regions directly; however, whether fore¬ 
ground residuals exist in the WMAP foreground-reduced maps is still an open 
question. Here, we use Pearson correlation coefficient (PCC) analysis with AS-7 
TeV cosmic ray (CR) data to probe possible foreground residuals in the WMAP 
nine-year data. The correlation results between the CR and foreground-contained 
maps (WMAP foreground-unreduced maps, WMAP template-based and MEM 
foreground models) suggest that: (1) CRs can trace foregrounds in the WMAP 
data; (2) at least some TeV CRs originate from the Milky Way; (3) foregrounds 
may be related to the existence of CR anisotropy (loss-cone and tail-in struc¬ 
tures); (4) there exist differences among different types of foregrounds in the 
declination range of < 15°. Then, we generate 10,000 mock CMB sky maps to 
describe the cosmic variance, which is used to measure the effect of the ffuctua- 
tions of all possible CMB maps to the correlations between CR and CMB maps. 
Finally, we do correlation analysis between the CR and WMAP foreground- 
reduced maps, and hnd that: (1) there are signihcant anticorrelations; and (2) 
the WMAP foreground-reduced maps are credible. However, the signihcant anti¬ 
correlations may be accidental, and the higher signal-to-noise ratio (SNR) Planck 
SMICA map cannot reject the hypothesis of accidental correlations. We therefore 
can only conclude that the foreground residuals exist with ~ 95% probability. 
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INTRODUCTION 


The Wilkinson Microwave Anisotropy Probe (hereafter WMAP) is designed to obtain 
a fnll sky map of the temperatnre anisotropy of the cosmic microwave backgronnd ra- 
diation (CMB) and has provided ns with qnite nsefnl information in cosmology research 


( Hinshaw et al.ll2013l) . The WMAP freqnency bands are 22, 33, 41, 61 and 94 GHz, named 
as K, Ka, Q, V, and W bands, and the corresponding angnlar resolntions are 0.88, 0.66, 
0.51, 0.35, and 0.22 deg, respectively. The WMAP team combines the CMB-dominated 
bands (Q, V and W bands) to obtain a signal-to-noise ratio (SNR) improved map; however, 
foregronnds as well as point sonrces also contribnte to the map at these selected bands. A 
linear combination of foregronnd templates is £t to the WMAP sky map at each freqnency 
band, and the £t is snbtracted to prodnce the rednced for egronnd map; the fo regronnd tem¬ 


plates are the WM AP K-Ka temperatnre difference map flHinshaw et al.ll2007r). the Ha map 


flFinkbeinerl l2003h and the FDS 94 GHz dnst emission map flFinkbeinerl Il999l) . Moreover, 


some techniqnes, inclnding internal linear combination (ILG), Maximnm Entropy Method 
(MEM), Markov Ghain Monte Garlo (MGMG) and are also implemented to characterize 
the natnre of the foregronnds (synchrotron, free-free, spinning dnst and thermal dnst). Then 
averaged foreground models (MEM, MGMGg, and Model 9) show that the total WMAP 



Ever since the release of the WMAP hrst-year data, the foreground residuals in the 
WMAP foreground-reduced data have received much attention. Joint multiwavelength cor¬ 
relation analysis is an important method fo r this purpose. Many kinds of data ha ve been used, 
such as the 10 and 15 GHz Tenerif e data (Ide Oliveira-Gosta et al.lll999L 120041) . the thermal 


ust and spinnin g dust templates flFinkbeinen 120041), the interstellar neu 


flVerschuuiil2007l) . the Sloan Digital Sky Surve y data flPadmanabhan et a 


20081) , the diffuse gamma-ray intensity maps (jWibig fc Wolfendale 
etc. 


2005 


ral hydrogen data 


2005 


Hirata et ah 


Liu fc Zhandl2006l) . 


Some basic components of the interstellar medium (ISM), which contain Galactic mag¬ 
netic held (GMF), Galactic cosmic ray (GGR) and the ordinary matter (such as gas and 
dust), have comparable pressures and are coupled together by electromagnetic forces. Both 
the dynamics of the ordinary matter and its spatial distribution at all scales are inhuenced by 
GGRs and GMFs; the GMFs and the propagation of GGRs provide efficient support against 
the self-gravity of ordinary matter. Gonversely, the inertia of the ordinary matter confines 
GMFs and, hence, GGRs, while its turbulent motions are responsible for the amplihcation of 
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GMFs and the acceleration of GCRs fjFerriefell200ll) . The WMAP foregronnd emissions are 
closely related to these components of the ISM: the rapid spir aling motions of GGR electron s 
abont magnetic field lines generate the synchrotron emissions flGinzburg fc Svrovatskiilll965h : 
the free-free emissions arise from ionized g as; and the spinning d nst an d thermal dnst emis¬ 
sions are obvionsly from interstellar dnst. IWibig fc Wolfendalel (120051) fonnd that there are 
foregronnd residnals in the WMAP first-ye ar data, which they interpreted as related to cos¬ 
mic ray (GR) intensity. ILin fc Zhang! (120061 1 also fonnd that there are residnals in the WMAP 
first-year foregronnd cleaned data. They correlated the WMAP first-year foregronnd-rednced 
data with the diffnse gamma-ray intensity maps from the Energetic Gamma-ray Experiment 
Telescope (EGRET), and fonnd that there are signihcant correlations at scales aronnd 15° 
in Q and W bands. They conclnded that the correlations are cansed by the foregronnd 
residnals in the WMAP data, which may be indnced by cosmic rays becanse at least some 
diffnse gamma-ray emissions are believed to be prodnced by the interactions of cosmic rays 
with gas and photon fields in the Milky Way. Hence, direct correlation analysis between 
the GR data and the WMAP data may be helpfnl to nnderstand the natnre of foregronnd 
residnals in the WMAP data. 

The WMAP team has released data several times; the latest and final prodncts are the 
WMAP nine-year data. The TeV GR intensity maps can be nsed to do correlation analysis 
with the WMAP nine-year data. Sever al TeV GR experimen ts have been c ondncted in recen t 


decades, incln ding Snper-Kamio kande (IGnillia 


ARG O-YBJ (IZhang et al.ll2010[) . EAS-TOP (lAglietta et al.N2009[) . IceGnbe (lAbbasi et ah 


1 et al.ll2nn7 ). MILA GRO (lAbdo et al.ll2nn9l). 


201 if) , etc. These experiments have published scientihc results, including GR anisotropy. 
The GR data we use here are from the large data sample of Tibet Air Shower Arrays (AS- 
7 ) experiment at 4, 6 , 12, 50 and 300 TeV, respectively; the arrays have been operated at 
Yangbajing (90.522° E, 30.102° N, 4,300 m above sea level) in Tibet, Ghina since 1990. The 
AS -7 experiment has an large effective area (up to 36,900 m^) and provides us enough GGR 
events (here ~ 37 billion s); here the experiment cannot distinguish gamma rays from the 
charged GR background (lAmenomori et al.ll2006l) . but cosmic rays dominate over gamma 
rays anyway in the Tibet air shower data. In our pa per the original GR rela tive intensity 
maps before smoothing are the same as in Figure 3 of lAmenomori et al.l (120061) . except that 
unity is subtracted from the intensity. There are two known st ructures of anisotropy on 


the GR map, i.e., the lo ss-cone structure and the tail-in structure (iNagashima fc Morilll976 


|P 


Nagashima et al.l 119981) . The loss-cone structure is the region with lower GR intensity and 
looks li ke a circle on the GR m ap, while the GR intensity has enhancement in the tail-in 


region. lAmenomori et al.l (120061 ) concluded that these anisotropic structures are not related 


to the heliospheric magnetic field but modulated by GMF, since they exist at energies of 
a few tens of TeV where the gyroradius is 1 AU. iGiacinti fc Sigll (120121) showed that 
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such structures arise from some structures in the tur bulent GMF within the CR scattering 
length (i.e., a few tens of parsecs from the Earth). iDrurv fc Aharonianl (120081) suggested 
that such structures could result from more distant sources if the magnetic mirrors allow for 
relatively scatter-free CR propagation from the source to the solar system. Hence, the CR 
anisotropy most likely originates from sonie loc al structures in close proximity to the solar 
system. In Figure 3 of lAmenomori et ahl (1200611 . there is a newly discovered structure, the 
Cygnus region. All these structures cha nge with the CR energy bands, and fade at higher CR 
energy bands ( Amenomori et aP 2006); for example, the CR maps at 300 TeV do not show 
any signihcant anisotropy. These anisotropic structures are important probes of foreground 
residuals in the WMAP data. 


Planck is the third generation of space-based CMB observatory operated by the Eu¬ 
ropean Space Agency, and has released its products in 2013. Compared with the WMAP, 
the Planck has wider frequency range (25-1000 GHz) and smal ler angular resolution (5-30 
arcmin utes), and provides higher SNR temperature sky map (jPlanck Collaboration et ah 


2013a|). The correlation analysis between the Planck data and the CR data may help to 


determine whether foreground residuals in the WMAP data exist. 


The paper is organized as follows: in Section |2] we describe the pre-processing of sky 
maps and the correlation analysis method. In Section |3] we prove that the CR data can be 
used to probe the WMAP foregrounds. In Section 0] we analyse the correlations between 
the CR data and the WMAP foreground-reduced data, and discuss the probabilities of 
foreground residuals’ existence. Finally, we summarize and make discussion in sectional 


2. CORRELATION ANALYSIS METHOD 

The correlation analysis between two maps needs a method to calculate correlation coef- 
hcients, and the Pearson correlation coefficient (PCC) method is commonly used. However, 
the sky maps detected by different instruments may be different in data formats and image 
qualities; for example, the WMAP maps at different frequency bands correspond to different 
angular resolutions. Before correlation analysis, we need to make maps with the same format 
and quality. 


2.1. Pre-Processing of Sky Maps 

We process all maps by three main steps, including coordinate transformation, map 
masking and smoothing. The AS-y CR sky maps sample data over square grids of size 
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with 2° in the declination range of about —10° to 70°, while other all-sky maps (WMAP, 
Planck and foregrounds) are stored with much smaller HEALPix grids. Hence, a practical 
option is that all maps should be transformed into the AS-y style map. In addition, the 
WMAP team uses masking templates to mask point sources and sky regions with uncertain 
fo regrounds; here the WMAP Extended Temperature Data Analysis Mask (i.e., the KQ75y9 


m 


Bennett et al.ll2013l) is chosen to cut as many uncertain regions as possible. After masking 


all maps, we find that the tail-in regions are partly masked except the declination range of 
less than 15° and the Cygnus regions are masked totally in the CR maps. Finally, all maps 
are smoothed to the same angular resolution (= 5°) by Guassian smoothing. The FWHM 
of the Gaussian smoothing is calculated by using 


FWHM = ^FWHM^ - FWHM^, (1) 

where FWHMq is the angular resolution to be smoothed to, FWHMj is the angular resolution 
of map i. 


2.2. Pearson Correlation Coefficient Method 


The PCG method is used to do correlation analysis. It is defined as 

n Ei XiVi - Ei Xi Ei Vi 


T = 
xy 




( 2 ) 


where x and y are the two data samples (e.g., the CR data and the WMAP data), n is the 
total number of elements in each sample, and rl.y is the PCG of x (CR) and y (WMAP). 
Masked pixels by KQ75y9 are not used in this correlation analysis. The values of PCCs 
range from -1 to 1, and there is no correlation when r\.y = 0. 


The errors related to rl.y are obtained by using 10,000 Monte Carlo simulations that 
considering the errors of x and y. r^y is transformed to z^y by the Fisher transformation 
equation: 


z 


I 

xy 


0.51n 


1 -I J..I 
' xy 

1 — 

xy 


( 3 ) 


Then z}^^y approximately follows a normal distribution with the standard deviation (i.e., 
the error of zL,). Hence, the error of the correlation coefficient ai can be calculated with 
the inverse Fisher transformation. 


By the above pre-processing and PCG method, we can obtain the curves of (between 
two sky maps) with declinations; here a correlation coefficient is determined independently 
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for each declination bin (~ 2°) becanse ab solnte calibration of the Tibet air shower detector 
efficiency in declination is not available flAmenomori et al.l l2006h . These cnrves actnally 
measnre the correlation of the two sky maps in the spatial distribntion, because the positive or 
negative values of r^y correspond to correlations or anticorrelations in the spatial distribution 
that the two maps have; here the correlation or anticorrelation increases with increasing 
absolute value of and the two maps have no correlation when rl.y = 0. Besides, if the sky 
maps for correlation analysis are isotropic, then the correlation coefficient r^y = 0; therefore, 
the changes of spatial distributions (e.g., anisotropic structures in the CR maps) are very 
important. There are two explanations for the observed correlation: two maps are all effected 
by the same physical processes, so the correlation is causal; or, if the two maps have no any 
physical relationship, the correlation is accidental. It is important to distinguish between 
these two possibilities. 


3. CORRELATION BETWEEN CR AND FOREGROUNDS 

In this section, we study if the CR data can be used as the probe of foreground residuals. 
We therefore need to do correlation analysis between the CR maps and some sky maps 
with foregrounds, i.e., the WMAP foreground-unreduced maps, the WMAP template-based 
foreground models and the MEM foreground models. 


3.1. Correlation betweeen CR and WMAP Data 

With the pre-processing described in Section [2Tl we can obtain the CR maps in different 
energy bands and the WMAP maps in different frequency bands, as shown in Figures [1] 
and [2j Then we use the CR maps to do correlation analysis with the WMAP foreground- 
unreduced maps and the total WMAP template-based foreground models with the PCC 
method described in Section and obtain the curves of PCCs and their errors as functions 
of declinations as shown in Figures [3] and 01 

The WMAP foreground-unreduced data contain the CMB and the total WMAP fore- 
gounds, and cover all hve frequency bands (K, Ka, Q, V, and W). From the curves of PCCs 
between the CR and WMAP foreground-unreduced data shown in Figures [3] and 01 we hnd 
that: a) These PCCs are positive in K band and fall to negative as frequencies increase till 
V and W bands. In the declination range of about 15° to 55°, where the CR maps only have 
the loss-cone structures, these correlations fall off very fast as frequencies increase; whereas 
in the declination range of about —10° to 15°, where the CR maps have both the tail-in and 
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loss-cone structures, they fall off more slowly as frequencies increase; b) These PCCs become 
weak for CR energy bands of 50 and 300 TeV. 

The total WMAP template-based foreground models only cover three frequency bands 
(Q, V, and W). From the curves of PCCs between the CR data and total WMAP template- 
based models, we find that: a) These PCCs are positive, and increase slightly with increasing 
WMAP frequencies till V and W bands; b) These PCCs become weak for CR energy bands 
of 300 TeV. 


The anisotropic structures (loss-cone and tail-in) in the CR maps fade with increasing 
energies; this explains why PCCs become weak at higher energies. The spectra of the 
declination-averaged correlation coefficients for the different CR energy bands (4, 6 and 12 
TeV) between the CR and WMAP foreground -unreduced data are similar to the spectra 
of total WMAP foregrounds (see Figure 22 in iBennett et ahl 120131) as shown in Figure El 
so the CR maps may be used to trace the foregrounds in the WMAP data. Moreover, the 
significant positive PCCs between the CR data and total WMAP template-based foreground 
models show that they have correlations of spatial distribution; these correlations are more 
likely caused by the same physical mechanisms rather than accident, because the TeV CR 
can be significantly influenced by the Galactic environment (e.g., magnetic held) and the 
WMAP foregrounds are mainly Galactic components. Hence, at least some TeV GR we 
observe should originate from the Milky Way. 


3.2. Correlation betweeen CR and MEM Data 


The total WMAP template-based foreground models only cover Q, V and W bands; 
however, the Maximum Entropy Method (MEM) foreground models (synchrotron, free-free, 
spinning dust and thermal dust) enable us to do correlation analysis between the CR and 
these foregrounds to cover all WMAP frequency bands. The MEM generates foregrounds by 
using the smoothed WMAP maps minus the ILC map as input and some fo reground maps as 
priors; the priors are the Haslam 408 MHz 


100 micron map fISchlegel et a 


thermal dust map ( Finkbeine: 


map 


the Ha map descr ibed iniFinkbeineii (20031) and iBennett et ahl (l2013l) for free-free, the SFD 


flHaslam et ah 


981 


19821) for synchrotron. 


19981) for spinning dust, and the FDS predicted 94 GHz 


199 


With the pre-processing described in Section 12.11 we 


can obtain the MEM foreground maps at K and Ka bands, as shown in Figured Then we 
use the CR maps to do correlation analysis with these MEM maps with the PCC method 
described in Section 12.21 and obtain the curves of PCCs and their errors as functions of 
declinations as shown in Figure [3 





































From the curves of PCCs between the CR and MEM maps, we hnd that: a) These PCCs 
are almost all positive except three foregrounds (synchrotron, spinning dust and thermal 
dust) in the declination range of —10° to 0°; b) These curves of PCCs at K band are 
identical with those at Ka band; c) The PCC curves of four types of foregrounds have 
relatively signihcant differences in the declination range of —10° to 15°; d) These PCCs 
become weak at 300 TeV. 

Similar to that in Section IXTl the PCCs become weak at high CR energy bands, due to 
the weak CR anisotropy at higher CR energies. In addition, the identical curves of PCCs at K 
and Ka bands suggest that each foreground’s spatial distribution remain unchanged. In fact 
the WMAP team generated maps of a type of MEM foreground at each frequency band by 
multiplying this type of foreground’s template by frequency dependent and pixel independent 
coefficient, i.e., there are the same spatial distribution for a type of MEM foreground map 
among different frequency bands; hence, the curves of PCCs at Q, V and W bands are 
identical with those at K or Ka band. Compared with the PCC curves of the declination 
range of > 15° (only cover the loss-cone structure), there are more signihcant differences in 
the declination range of < 15° (cover the loss-cone and tail-in structures) among four types 
of foregrounds; the loss-cone and tail-in features may have some physical differences, or the 
loss-cone features are different between these two declination ranges. 

In Section 13.11 the PCC curves of the total WMAP template-based foreground models 
with the CR maps are slightly different among Q, V and W bands. Actually, each type of 
WMAP template-based foreground at all frequency bands also have identical spatial distri¬ 
bution; however, the proportion of different type of foreground in the total foreground model 
changes with increasing frequency, so the spatial distribution of the total foreground model 
changes as well. The PCCs of the foreground models (MEM and WMAP template-based) 
with the CR maps are almost all signihcantly positive; these correlations are more likely 
causal and can be used to understand the WMAP foregrounds because the TeV CR and 
the WMAP foregrounds are all related to Galactic ISM. Hence, CRs may be a good tracer 
of foregrounds at all WMAP frequencies. Combined with the PCC curves of the WMAP 
foreground-unreduced data with the CR maps in Section 13.11 we can conclude that the total 
foregrounds in the WMAP foreground-unreduced data dominate at K and Ka bands while 
the CMB dominate at Q, V and W bands; hence, the spectra of the declination-averaged 
correlation coefficients agree with the spectra of total foregrounds. Moreover, the results of 
this section also suggest that the values of PCCs between the CMB and CR maps may be 
negative. 
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4. FOREGROUND RESIDUALS IN WMAP DATA 

Section |3] shows that the CR data are useful for probing the WMAP foregrounds, so 
we can do correlation analysis with the WMAP foreground-reduced data. However, possible 
accidental correlations between the CR and CMB sky maps have an effect on hnding out the 
foreground residuals, due to cosmic variance. Hence, we need theoretical CMB sky maps to 
estimate possible impacts of accidental correlations, and the latest Planck sky map is also 
helpful. 


4.1. Cosmic Variance 


The CMB sky map detected by the WMAP is one sample of the population including 
many possible CMB sky maps (e.g., observed from other galaxies) with identical cosmological 
parameters; therefore, the cosmic variance is used to measure the effect of the fluctuations 
of all possible CMB maps to the correlations between CR data and CMB sky maps. Hence, 
we need large amounts of mock CMB sky maps to reflect intrinsic fluctuation of CMB in 
spatial distribution and obtain the cosmic variance. 


SYNFAST, i.e, a program in HEALPix software, is used to create 10,000 mock CMB sky 
maps (smoothed to 5°) with the same angular power spectra generated by CAMB tool; here 
we use the cosm ological parameters p rovided by the WMAP nine-year cosmological results 
(see Table 17 in iBennett et al.ll2013l) . Every mock CMB sky map has a random spatial 
distribution of temperature, and then 10,000 correlation coefficients (r^^., calculated with 


xy 

Equation ([2])) are calculated between the CR data and 10,000 mock CMB sky maps at each 
declination. Finally, can be transformed to with Equation (|3]), and then 10,000 
each declination approximately follow a normal distribution with the mean < 2 :: 


xy 

> and the 


standard deviation Therefore, the mean < >, the one-sigma and the two-sigma 

can be derived from < , and 2 by the inverse Fisher transformation, 

respectively; here the cosmic variance can be represented by (7^1^, while represent the 
measurement variance of an instrument (e.g., the detectors in the WMAP observatory). 


If there are only accidental correlations, < > approaches to zero with increasing 

sample numbers; otherwise, the CR and CMB maps have causal correlations. In addition, if 
we want to ensure that the observed CMB sky map is credible, the measurement variance 
must be much smaller than the cosmic variance (i.e., aj.1^ 3> Figure [8] shows the 

statistics results of 10,000 sets of PCCs between the CR and mock CMB maps: the mean 
< ^^xy >) lines and 2-a (ct”^, 2 <x) lines. We hnd that all mean values (i.e., < >) 

are close to zero; therefore, we conclude that the CR and CMB maps only have accidental 
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correlations, as expected. 


4.2. Correlation betweeen CR and WMAP Foreground-reduced Data 

The WMAP foreground-reduced data only cover Q, V and W bands, and are used to 
do correlation analysis with the CR maps with the PCC method described in Section 12.21 
The curves of PCCs and their errors with declinations are shown in Figure [9l here l-cx lines 
and 2-a lines are drawn in hgures. 

For the curves of PCCs between the CR and WMAP foreground-reduced data, we hnd 
that: a) b) These PCCs are negative and signihcant (Ir^yl > from 4 to 12 

TeV in the declination range of about 15° to 35°, and become weak for CR energy bands of 
50 and 300 TeV; c) These PCCs have almost no change from Q to W bands at the same CR 
energy bands. 

Therefore, we can conclude that: (1) the measurement variance is much smaller than 
the cosmic variance, as described in Section 14.11 and (2) there are anticorrelations between 
the CR and WMAP foreground-reduced maps. 


4.3. Correlation betweeen CR and Planck Data 


The Planck space mission has released higher SNR maps (jpianck Collaboration et~ah 


2013alJbl) . so that we can compare the WMAP ILC sky map with the Planck SMICA sky 


map. These are CMB maps formed by linearly combining data from all WMAP bands or 
all Planck channels with weighting such that foreground emission is approximately canceled. 
The causal correlations with CR data caused by foreground residuals may result in different 
PCCs of these two maps. The WMAP ILC and Planck SMICA maps are processed with the 
pre-processing described in Section |2T] as shown in Figure flOl and the curves of PCCs and 
their errors with declinations are shown in Figure fTTl 


In Figure dH we hnd that: a) b) These PCCs are negative and signihcant 

(kxyl > ^rly) from 4 to 12 TeV in the declination range of about 15° to 35°, and become 
weak for CR energy bands of 50 and 300 TeV. Therefore, we can conclude that: (1) the 
measurement variance is much smaller than the cosmic variance, as described in Section ITTl 
and (2) there are anticorrelations between the CR and WMAP ILC (or Planck SMICA) 
maps. 


The PCC curves of both WMAP ILC and Planck SMICA maps have no signihcant 
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difference, and the PCCs with the CR data of the WMAP ILC map are stronger slightly 
than those of the Planck SMICA map in the declination range of abont 15° to 35°. Hence, we 
cannot conclnde that the Planck SMICA map have more serions problems abont foregronnd 
residnals than the WMAP ILC map. 


4.4. Probability of Accidental Correlation 


The rednced foregronnd sky maps (WMAP foregronnd-rednced, WMAP ILC and Planck 
SMICA maps) and the CR maps have signihcant {\rl.y\ > (rlly) anticorrelations from 4 to 12 
TeV in the declination range of abont 15° to 35°; here the loss-cone strnctnres are located in 
this declination range. However, we cannot reject the hypothesis that these anticorrelations 
are accidental, before we make qnantitative analysis to calcnlate probabilities that these 
anticorrelations are accidental. 


There are 10,000 mock CMB sky maps in Section ITTl and we can find ont the nnmber 
that the correlations of the mock CMB and CR maps are stronger than those of the rednced 
foreground and CR maps, i.e., the absolute value of the sum of weighted PCCs between 
the mock CMB and CR maps | the absolute value of the sum of weighted 

PCCs between the reduced foreground and CR maps satisfy: — 

\J2''^xy^xy\] here the weight W^y (oc cos(Dec.)) is used to account for the variation of 
the solid angle of the sampled sky with declinations. Then the probability of accidental 
correlation (AC) can be calculated by 


Pac 


Nac 

10 , 000 ’ 


(4) 


where Nac is the number of maps with | ^ i^^yW^yl > \ '^ly^xyl among 10,000 mock CMB 
sky maps in Section 14.11 


We need to choose a reasonable declination range before calculating Pac- Full available 
declination range can be used to describe accidental correlation in the declination range of 
about —10° to 70°; however, narrower declination range is more reasonable when foreground 
residuals only dominate at local ranges. Therefore, the significant (|r^y| > anticorre¬ 
lations in the declination range of about 15° to 35° are worth special attention. Table [1] 
lists Pac of WMAP (Q to W band and ILC) and Planck SMICA maps from 4 to 12 TeV 
in three different declination ranges: ‘A’ for full available declination range, ‘B’ for —10° 
to 35° declination range (cover the loss-cone and tail-in structures) and ‘C’ for 15° to 35° 
declination range (only cover the loss-cone structure); here the probabilities of 50 and 300 
TeV CR energy bands are not calculated, because the CR maps at these energy bands are 
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nearly isotropic. We can conclude that: (1) the foreground residuals exist in the WMAP 
foreground-reduced maps with ~ 95% probability at scales of 20°, but accidental correlations 
with ~ 5% probability cannot be excluded; (2) compared with the WMAP ILC map, the 
Planck SMICA map has lower foreground residual levels because of the larger Pag; (3) the 
data in the declination range of < 15° have great impact on the probabilities, because ‘B’ 
range’s maximum Pag difference (23%) is much larger than ‘A’ range’s (7%) or ‘C’ range’s 
(10%), i.e., ‘B’ range’s Pag has the biggest change with the WMAP frequency or CR energy 
bands; (4) the values of Pag decrease with increasing WMAP frequency or CR energy bands. 


5. SUMMARY AND DISCUSSION 

We have first performed correlation analysis of the AS-y CR data with the foreground- 
contained sky maps (i.e., the WMAP foreground-unreduced maps, the WMAP template- 
based foreground models and the MEM foreground models) to conhrm that the CR can be 
used to probe foregrounds in the WMAP data. Next, we generate 10,000 mock CMB sky 
maps to describe the cosmic variance, and then analyzed the signihcant anticorrelations of 
the CR data with the WMAP foreground-reduced maps. Moreover, the latest Planck SMICA 
map is used to compare with the WMAP ILC map. In order to understand whether the 
foreground residuals exist, we obtain the probabilities of accidental or causal correlations. 
The main conclusions of this paper are as following. 

1) The spectra of the declination-averaged correlation coefficients for the different CR 
energy bands (4, 6 and 12 TeV) between the CR and WMAP foreground-unreduced data is 
similar to the spectra of total WMAP foregrounds. Hence, the CR data can be used to trace 
the WMAP foregrounds. 

2) The correlations of the CR data with the WMAP template-based foreground models 
and the MEM foreground models (synchrotron, free-free, spinning dust and thermal dust) 
suggest the Galactic origin of at least some TeV CRs. 

3) The correlations of the WMAP foregrounds with the anisotropic structures (loss-cone 
and tail-in structures) in the CR maps show that the existence of these structures may be 
related to some foreground emissions (e.g., synchrotron emissions). 

4) The correlations of the CR data with four types of foregrounds have relatively sig¬ 
nihcant differences in the declination range of < 15°, and the maximum Pag difference of 
—10° to 35° declination range (23%) is much larger than the maximum Pag difference of 15° 
to 35° declination range (10%); hence, the CR data in the declination range of < 15° have 
great impact on the correlation results. 



5) The cosmic variances are far greater than the measurement variances, so the WMAP 
foreground-reduced maps are credible for making correlation analysis with the CR data. 

6) The foreground residuals exist in the WMAP foreground-reduced maps with high 
probabilities ~ 95% at scales of 20°, but accidental correlations ~ 5% cannot be excluded. 

7) The Planck SMICA map has lower foreground residual levels than the WMAP ILC 
map, because the accidental correlation probabilities of the former are greater than those of 
the latter. 

Three main components, including CMB, CR and foregrounds, are studied in this work. 
Since the total WMAP template-based foreground models cannot contain all foreground 
emissions, so some foregrounds may still exist in the WMAP foreground-reduced map, i.e, 
the hnal WMAP CMB map may not be clean. Hence, the CR data are used to hnd out 
the foreground residuals, because the CR data we use have two advantages: (1) the TeV 
CR have close relationship with the foregrounds; and (2) the TeV CR have no any physical 
relationship with the CMB from the early universe at redshift z ~ 1100. However, the 
accidental correlations due to cosmic variance between the CR and CMB maps (i.e., only 
the correlations of spatial distributions, not a physical relationship) have an effect on probing 
the foregrounds; nevertheless we can calculate the probabilities of the accidental (or causal) 
correlations. 

The anisotropic structures (loss-cone and tail-in structures) play an important role in 
the correlation analysis. For both the total WMAP template-based foreground models and 
the MEM foreground models (synchrotron, free-free, spinning dust and thermal dust), the 
correlations with the CR maps are related to the loss-cone structure and the tail-in structure, 
i.e., the foreground maps have some similar structures. These structures may result from 
(1) the propagation of the cosmic rays in some regions of the Milky Way changes the spatial 
distribution of foreground emissions; and (2) the spatial distribution of the CRs and the 
foregrounds are all affected by some common physical mechanisms. However, cosmic rays are 
not directly involved in free-free, spinning dust and thermal dust emission mechanisms, and 
even the synchrotron emission associated with CR anisotropy would be negligible compared 
to the synchrotron emission from the general ISM; here CR anisotropy is most likely the local 
features (near the solar system), so the associated synchrotron emission originate over a short 
path length rather than entire path length through the Galaxy. One possible explanation 
is that: there exist many CR producing anisotropic structures throughout the Milky Way, 
those distant structures produce a roughly isotropic CR map while some local structures lead 
to anisotropy in the CR map. Then these local structures, which contain gas and dust, also 
generate all kinds of microwave foregrounds correlated with CR anisotropy. Moreover, each 
anisotropic structure has its own characteristics, including motions and components (gas. 
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dust, magnetic filed, supernova remnant, etc.), and then different types of structures may 
have different effect (inhibition or promotion) on CR propagation. Hence, the distinction 
between the tail-in and loss-cone structures may be due to different types of anisotropic 
structures. 

Whether the foreground residuals exist in the WMAP foreground-reduced map has not 
been solved completely, because we cannot exclude the accidental correlations with high 
conhdence. Hence, further works are needed: a) better foreground templates and higher 
SNR CMB sky map are helpful; b) the formation mechanism of the anisotropic structures 
in the CR map is still an open question, and further studies about the loss-cone and tail-in 
structures may be useful to understand our correlation results; c) more available absolute 
calibration of CR detectors and more CR data over larger sky areas can help us to better 
probe the foreground residuals. 
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Table 1: Probabilities of accidental correlations between the WMAP Foreground-reduced 
Maps and the CR Maps 



Q band 

V band 

W band 

ILC 

SMICA 


(%) 

(%) 

(%) 

(%) 

(%) 

A: full available declination range 

4 TeV 

ST.rtri 

36.911^ 

32.211? 

41.811® 


6 TeV 


35.11®;® 

30.711? 

41.8l^°2® 


12 TeV 

36.8l^°4' 

36.11^1' 

32.2l|? 

43.51?°8® 


B: - 

.0° to 35° declination 

range 


4 TeV 

44.71^;^ 

43.0l|e^ 

36.011® 

50.41IJ 

59.411? 

6 TeV 

33.611® 

31.81®:® 

26.011^ 

40.21” 

44.611^ 

12 TeV 

29.111® 

27.6l|o" 

22.411? 

qc 1 +6.9 
'^'-'•^--6.0 

43.811? 


C: 15° to 35° declination 

range 


4 TeV 

14.811^ 

14.111^ 

12.111® 

16.211® 

23.711® 

6 TeV 

8.01};® 

7.61};^ 

6.311® 

9.111^ 

12.91?;? 

12 TeV 

6.2l}l 

5.81};® 

4.611^ 

y 0+1-9 

1 1 7+2-5 
'-1.8 
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AS- Y Intensity Map 
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Fig. 1.— Smoothed AS -7 cosmic ray (CR) intensity maps and WMAP nine-year foregronnd- 
unreduced temperatnre flnctnation maps in eqnatorial coordinates, sampled over a sqnare 
grid of size width 2°; all maps are smoothed to 5°. The circled regions labeled by I and 
II in the left top three panels are the tail-in regions and the loss-cone regions, respectively. 
The right hve panels are smoothed WMAP nine-year foregronnd-nnrednced temperatnre 
flnctnation maps, which are masked by the KQ75y9 and the point sonrce catalog masks, in 
nnits of mK; these masks are shown in dark grey in the WMAP foregronnd-nnrednced maps 
and contonrs of their bonndaries are added to the CR maps. Dnring the transformation 
from Galactic coordinates to equatorial coordinates, a 2 ° x 2 ° grid will be set to zero and 
not be used in the correlation analysis if more than 50% of HEALPix pixels in the grid are 
masked. 
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WMAP foreground-reduced map WMAP template-based foreground model 
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Fig. 2.— Smoothed WMAP nine-year foreground-reduced temperature fluctuation maps 
and WMAP template-based foreground models at Q, V, W bands, in units of mK. The 
smoothing and mask procedures are the same as Figured! 
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Fig. 3.— Correlation coefficients (green triangles) of the AS -7 intensity maps (the left hve 
panels in Fignre [T]) with the WMAP nine-year foregronnd-unrednced data at K and Ka 
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Fig. 4.— Similar to Figure |3l correlation coefficients (green triangles and blue squares) of the 
AS -7 intensity maps (the left hve panels in Figured]) with the WMAP nine-year foreground- 
unreduced data at Q ,V and W bands (the right bottom three panels in Figured]) and the 
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Fig. 5.— Spectra of the declination-averaged correlation coefficients between the WMAP 
foreground-unrednced data and the CR maps at 4, 6 and 12 TeV (green triangles, red circles 
and blue squares), and spectra of total foreground models (black solid curves) outside of the 
KQ75y9 mask. 
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Fig. 6.— Maximum Entropy Method (MEM) maps (synchrotron, free-free, spinning dust, 
and thermal dust foregrounds) at K and Ka bands. These maps are smoothed to an angular 
resolution of 5°. 
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Fig. 7.— Similar to Figure [31 correlation coefficients (red circles, green triangles, blue squares 
and yellow diamonds) of the AS-y intensity maps (the left hve panels in Figure [1]) with the 
synchrotron, free-free, spinning dust, and thermal dust MEM maps at K and Ka bands (the 
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Fig. 8.— The statistical results of 10,000 sets of correlation coefficients between the CR 
maps and the mock CMB maps: the mean values (green triangles), l-cx cosmic variance lines 
(dotted lines) and 2-cr cosmic variance lines (dot-dashed lines). 
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Fig. 9.— Correlation coefficients (red circles) of the AS-7 intensity maps (the left five panels 
in Figure [I]) with the WMAP nine-year foreground-reduced data at Q, V and W bands 
(the left three panels in Figure [2]) as functions of declinations. We use l-cr cosmic variance 
lines (dotted lines) to divide non-significant correlations from significant correlations and 
2-cr cosmic variance lines (dot-dashed lines) to divide significant correlations from highly 
significant correlations. 
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Fig. 10.— The WMAP ILC sky map (reduced-galaxy map produced via linear combination 
technique) and the Planck SMICA sky map (reconstructed CMB map as a linear combina¬ 
tion). These two maps are smoothed to an angular resolution of 5°. 
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Fig. 11.— Similar to Figure [9l correlation coefficients (blue squares and red triangles) of the 
AS -7 intensity maps (the left hve panels in Figure [1]) with the WMAP ILC sky map (the 
left panel in Figure [m and the Planck SMICA sky map (the right panel in Figure fTOll . 

















